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Three different special quasirandom structures (SQS) of the substitutional hep Ai-^Bj: binary 
random solutions [x — 0.25, 0.5, and 0.75) are presented. These structures are able to mimic the 
most important pair and multi-site correlation functions corresponding to perfectly random hep 
solutions at those compositions. Due to the relatively small size of the generated structures, they 
can be used to calculate the properties of random hep alloys via first-principles methods. The 
structures are relaxed in order to find their lowest energy configurations at each composition. In 
some cases, it was found that full relaxation resulted in complete loss of their parental symmetry 
as hep so geometry optimizations in which no local relaxations are allowed were also performed. 
In general, the first-principles results for the seven binary systems (Cd-Mg, Mg-Zr, Al-Mg, Mo-Ru, 
Hf-Ti, Hf-Zr, and Ti-Zr) show good agreement with both formation enthalpy and lattice parameters 
measurements from experiments. It is concluded that the SQS's presented in this work can be widely 
used to study the behavior of random hep solutions. 
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I. INTRODUCTION 



Thermodynamic modeling using the CALPHAD 
(CALculation of PHAse Diagrams) methodii^ attempts 
to describe the Gibbs energy of a system through empiri- 
cal models whose parameters are fitted using experimen- 
tal information. These descriptions allow the extrapola- 
tion of a system's thermodynamic properties to regions in 
the composition-temperature space that have not/cannot 
be accessed through experiments. These empirical mod- 
els, however, are as good as the data used to fit them and 
are therefore limited by the availability of accurate exper- 
imental data. This limitation can be overcome by using 
theoretical calculations based on first-principles methods, 
which are capable of predicting the physical properties 
of phases with no experimental input"^. Unfortunately, 
despite their predictive nature, these methods are not 
yet able to calculate the thermochemistry of materials — 
especially multi-component, multi-phase systems — with 
the precision required in industry. 

A natural way of improving the predictive capabilities 
of empirical models while maintaining their applicabil- 
ity to practical problems is by combining first-principles 
and CALPHAD techniques. Thanks to efhcient schemes 
for implementing Density Functional Theory (DFT)^, 
the almost-routine use of first-principles results within 
the CALPHAD methodology has become a reality. In 
this hybrid approach, the energetics obtained through 
electronic structure calculations are used as input data 
within the CALPHAD formalism to obtain the parame- 
ters that describe the Gibbs energy of the system^. 



The first-principles electronic structure calculations 
of perfectly ordered periodic structures are relatively 
straightforward since they usually rely on the use of pe- 
riodic boundary conditions. Problems arise, however, 
when attempting to use these methods to study the ther- 
mochemical properties of random solid solutions since an 
approximation must be made in order to simulate a ran- 
dom atomic configuration through a periodic structure. 
The usual approaches that have been used in the past 
can be summarized as follows: 

i) The most direct approach is the supercell method. 
In this case, the sites of the supercell can be ran- 
domly occupied by either A or B atoms to yield 
the desired Ai-xBx composition. In order to repro- 
duce the statistics corresponding to a random alloy, 
such supercells must necessarily be very large. This 
approach is, therefore, computationally prohibitive 
when the size of the supercell is on the order of hun- 
dreds of atoms. 

ii) Another technique, the Coherent Potential Approx- 
imation (CPA)^ method, is a single-site approxima- 
tion that models the random alloy as an ordered lat- 
tice of effective atoms. These are constructed from 
the criterion that the average scattering of electrons 
off the alloy components should vanish?.. In this 
method, local relaxations are not considered explic- 
itly and the effects of alloying on the distribution 
of local environments cannot be taken into account. 
Local relaxations have been shown to significantly 
affect the properties of random solutions^, especially 
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when the constituent atoms vary greatly in size and, 
therefore, their omission constitutes a major draw- 
back. Although the local relaxation energy can be 
taken into account^, these corrections rely on cluster 
expansions of the relaxation energy of ordered struc- 
tures and the distribution of local environments is 
not explicitly considered. Additionally, such correc- 
tions are system specific. 

iii) A third option is to apply the Cluster Expansion 
approach^. In this case, a generalized Ising model 
is used and the spin variables can be related to the 
occupation of either atom A or B in the parent lat- 
tice. In order to obtain an expression for the config- 
urational energy of the solid phase, the energies of 
multiple configurations (typically in the order of a 
few dozens) based on the parent lattice must be cal- 
culated to obtain the parameters that describe the 
energy of any given Ai^^Bx composition. This ap- 
proach typically relies on the calculation of the en- 
ergies of a few dozen ordered structures. 

In the techniques outlined above, there are serious 
limitations in terms of either the computing power re- 
quired (supercells, cluster expansion) or the ability to 
accurately represent the local environments of random 
solutions (CPA) . Ideally, one would like to be able to ac- 
curately calculate the thermodynamic and physical prop- 
erties of a random solution with as small a supercell as 
possible so that accurate first-principles methods can be 
applied. This has become possible thanks to the devel- 
opment of Special Quasirandom Structures (SQS). 

The concept of SQS was first developed by Zunger 
et al.— to mimic random solutions without generating a 
large supercell or using many configurations. The basic 
idea consists of creating a small — 4-48 atoms — periodic 
structure with the target composition that best satisfies 
the pair and multi-site correlation functions correspond- 
ing to a random alloy, up to a certain coordination shell. 
Upon relaxation, the atoms in the structure are displaced 
away from their equilibrium positions, creating a distri- 
bution of local environments that can be considered to 
be representative of a random solution, at least up to the 
first few coordination shells. 

Provided the interatomic electronic interactions in 
a given system are relatively short-range, the first- 
principles calculations of the properties of these designed 
supercells can be expected to yield sensible results, espe- 
cially when calculating properties that are mostly depen- 
dent on the local atomic arrangements, such as enthalpy 
of mixing, charge transfer, local relaxations, and so forth. 
It is important to stress that the approach fails whenever 
a property depends on long-range interactions. 

The SQS's for fcc-based alloys and bcc alloys have been 
generated by Wei et al.— and Jiang et al.— , respectively. 
However, to the best knowledge of these authors, there 
has been no investigation on the application of the SQS 
approach to the study of hep substitutional random solu- 
tions. In the present work, we propose two SQS's capable 



of mimicking hep random alloys at 25, 50 and 75 at.%. 
The paper is organized as follows: 

The proposed SQS are characterized in terms of their 
ability to reproduce the pair and multi-site correlation 
functions of a truly random hep solution. Subsequently, 
the structures are tested in terms of their ability to re- 
produce, via first-principles calculations, the properties 
of certain selected stable or metastable binary hep so- 
lutions, namely Cd-Mg, Mg-Zr, Al-Mg, Mo-Ru, Hf-Ti, 
Hf-Zr and Ti-Zr. To further analyze the relaxation be- 
havior of the structures, the distribution of first nearest 
bond lengths as well as the radial distribution for the first 
few coordination shells is presented. Finally, for each of 
the selected binaries, the calculated and available exper- 
imental lattice parameters and enthalpy of mixing are 
compared. Results from other techniques are also pre- 
sented where available in order to further corroborate 
the present calculations. 

II. GENERATION OF SPECIAL 
QUASIRANDOM STRUCTURES 

In order to characterize the statistics of a given atomic 
arrangement, one can use its correlation functioni^. 
Within the context of lattice algebra, we can assign a 
"spin value", cr = ±1, to each of the sites of the configu- 
ration, depending on whether the site is occupied by A or 
B-type atoms. Furthermore, all the sites can be grouped 
in figures, f {k,m), of k vertices, where A:=l,2,3,- ■ • re- 
sponds to a shape, point, pair, and triplet,. . . respec- 
tively, spanning a maximum distance of m, where m = 
1,2,3, ••• is the first, second, and third-nearest neigh- 
bors, and so forth. The correlation functions, Tik,m-i are 
the averages of the products of site occupations (±1 for 
binary alloys and ±1,0 for ternary alloys) of figure k at 
a distance m and are useful in describing the atomic dis- 
tribution. The optimum SQS for a given composition is 
the one that best satisfies the condition: 

(nfc,™)5Q5 - (^k.n)B. (1) 

where (Jlk.rnjR is the correlation function of a random al- 
loy, which is simply by {2x — 1)^ in the Ai-xBx substitu- 
tional binary alloy, where x is the composition. We con- 
sidered SQS's of two different compositions, i.e. x — 0.5 
and 0.75. 

Unlike cubic structures, the order of a given configu- 
ration in the hep lattices relative to a given lattice site 
may be altered with the variation of c/a ratio. However, 
these new arrangements will not cause any change in the 
correlation functions, since one can thus use any c/a ra- 
tio to generate the hep SQS's. As a matter of simplicity, 
the ideal c/a ratio was considered in order to generate 
SQS's. 

In the present work, we used the Alloy Theoretic Au- 
tomation Toolkit (ATAT)'^ to generate special quasiran- 
dom structures for the hep structure of 8 and 16 sites. 
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The schematic diagrams of the created special quasiran- 
dom structure with 16 atoms are shown in FIG. [1] and 
the corresponding lattice vectors and atomic positions 
are hsted in TABLE [J 

The correlation functions of the generated 8 and 16- 
atom SQS's were investigated to verify that they satisfied 
at least the short-range statistics of an hep random solu- 
tion. As is shown in TABLE |TT1 the 16-atom structures 
satisfy the pair correlation functions of random alloys 
up to the fifth and third nearest neighbor for the 50 at.% 
and the 75 at.% compositions, respectively. On the other 
hand, TABLE HI] shows that the SQS-8 for 75 at.% could 
not satisfy the random correlation function even for the 
first-nearest neighbor pair. Thus, SQS's with 16 atoms 
are capable of mimicking a random hep configuration be- 
yond the first coordination shell. 

It is important to note that in TABLE |TT1 and contrary 
to what is observed in the SQS for cubic structures, some 
figures have more than one crystallographically inequiv- 
alent figure at the same distance. For example, in the 
case of hep lattices with the ideal c/a ratio, two pairs 
may have the same interatomic distance and yet be crys- 
tallographically inequivalent. In this case, despite the 
fact that the two pairs (0,0,0) and (a, 0,0); (0,0,0) and 
(i,|,i), have the same inter-atomic distance, a, they do 
not share the same symmetry operations. This degener- 
acy is broken when the c/a ratio deviates from its ideal 
value. 

For the sake of efficiency, the initial lattice parameters 
of the SQS's were determined from Vegard's law. By do- 
ing so, the c/a ratio was no longer ideal. Afterwards, we 
checked the correlation functions of the new structures 
and found that they remained the same as long as the 
corresponding figures were identical. 

The maximum range over which the correlation func- 
tion of an SQS mimics that of a random alloy can be 
increased by increasing the supercell size. As the size 
of the SQS increases, the probability of finding config- 
urations that mimic random alloys over a wider coor- 
dination range increases accordingly. The search algo- 
rithm used in this work consists of enumerating every 
possible supercell of a given volume and for each su- 
percell, enumerating every possible atomic configuration. 
For each configuration, the correlation functions of dif- 
ferent figures, i.e. points, pairs, and triplets, are cal- 
culated. To save time, the calculation of the correla- 
tions is stopped as soon as one of them does not match 
the random state value. This algorithm becomes pro- 
hibitively expensive very rapidly. The generation of a 
larger SQS could be accomplished by using a Monte- 
Carlo-like scheme (e.g., Abrikosov et al. 1^), but this is 
beyond the scope of present work. In fact, the authors 
could generate a 32-atom SQS's, and the average total 
energy difference between 16-atom SQS's and 32-atom 
SQS's in the Cd-Mg system was around 2 meV per atom. 
The authors maintain a focus on 16-atom SQS, because 
this size represents a good compromise between accuracy 
and the computational requirements associated with the 



TABLE I: Structural descriptions of the SQS-A'^ structures 
for the binary hep solid solution. Lattice vectors and atomic 
positions are given in fractional coordinates of hep lattice. 
Atomic positions are given for the ideal, unrelaxed hep sites. 
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necessary first-principles calculations. 

It is also important to note that finding a good hep 
SQS is more difficult than finding an SQS of cubic struc- 
tures with the same range of matching correlations due to 
the fact that, for a given range of correlations, there are 
more symmetrically distinct correlations to match. Addi- 
tionally, the lower symmetry of the hep structure implies 
that there are also many more candidate configurations 
to search through in order to find a satisfactory SQS. 
Thus, the number of distinct supercells is larger and the 
number of symmetrically distinct atomic configurations 



4 




(a) SQS-16 for 2^=0.5 (b) SQS-16 for x=0.75 



FIG. 1: Crystal structures of the Ai-^Bx binary hep SQS-16 structures in their ideal, unrelaxed forms. All the atoms are at 
the ideal hep sites, even though both structures have the space group, PI. 



TABLE II: Pair and multi-site correlation functions of SQS-A*' structures when the c/a ratio is ideal. The number in the square 
bracket next to Ilk,m is the number of equivalent figures at the same distance in the structure, the so-called degeneracy factor. 
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is larger, in comparison to fee or bee lattiees. 

In order to verify the proposed 16-atom SQS's are ad- 
equate for the simulation of hep random solutions, the 
authors ealeulated other SQS's at 75 at.% which have 
random-like pair correlations up to the third nearest- 
neighbor but that have slightly different correlations for 
the fourth nearest-neighbor. The pair correlation func- 
tion at 75 at.% of a truly random solution would be 
(2 X 0.75 - 1)2 0.25 and therefore the four SQS's in 
TABLE mil are worse than the one used in the present 
work. These structures were applied to the Cd 25 at.%- 



Mg 75 at.% system and, as can be seen in TABLE IIIIl 
the associated energy differences are negligible. This is 
due to the fact that the energetics of this system are 
dominated by short-range interactions. Thus, as long 
as the most important pair correlations (up to the third 
nearest-neighbors in hep structure with ideal c/a ratio) 
are satisfied, the SQS's can successfully be applied to ac- 
quire properties of random solutions in which short-range 
interactions dominate. 
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TABLE III: Pair correlation functions up to the fifth and 
the calculated total energies of other 16 atoms sqs's for 
Cdo.25Mgo.75 are enumerated to be compared with the one 
used in this work (SQS-16). The total energies are given in 
the unit, eV/atom. 
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III. FIRST-PRINCIPLES METHODOLOGY 

The selected hep SQS-16 structures were used as ge- 
ometrical input for the first-principles calculations. The 
Vienna Ah initio Simulation Package (VASP)i^ was used 
to perform the Density Functional Theory (DFT) elec- 
tronic structure calculations. The projector augmented 
wave (PAW) methodic was chosen and the general gra- 
dient approximation (GGA)^^ was used to take into 
account exchange and correlation contributions to the 
hamiltonian of the ion-electron system. A constant en- 
ergy cutoff of 350 eV was used for all the structures, 
with 5,000 fc-points per reciprocal atom based on the 
Monkhorst-Pack scheme for the Brillouin-zone integra- 
tions. The k-point meshes were centered at the F point. 
The convergence criterion for the calculations was 10 
meV with respect to the 16 atoms. Spin-polarization was 
not taken into account. The generated SQS's were either 
fully relaxed, or relaxed without allowing local ion relax- 
ations, i.e., only volume and c/a ratio were optimized. 
As will be seen below, the full relaxation caused some of 
the SQS's to lose the original hep symmetry. 



IV. RESULTS AND DISCUSSIONS 

A. Analysis of Relsixed Structures 

The symmetry of the resulting SQS was checked using 
the PLATON^^ code before and after the relaxations. 
Both SQS's have the lowest symmetry of PI, although 
all the atoms are sitting on the lattice sites of hep. The 
procedure was verified by checking the symmetries of the 
generated unrelaxed SQS. Once all the sites in the SQS 
were substituted with one single atomic species, PLA- 
TON identified SQS's as perfect hep structures. All the 
atoms of the initial structures are on their exact hep lat- 
tice sites. However, upon relaxation the atoms may be 



displaced from these ideal positions. According to the 
definition of an hep random solution, all the atoms, in 
this case two different type of atoms, should be at the 
hep lattice points — within a certain tolerance — even 
after the structure has been fully relaxed. The default 
tolerance of detecting the symmetry of the relaxed struc- 
tures allowed the atoms to deviate from their original 
lattice sites by up to 20%. 

In principle, relaxations should be performed with re- 
spect to the degrees of freedom consistent with the initial 
symmetry of any given configuration. In the particular 
case of the hep SQS's, local relaxations may in some cases 
be so large that the character of the underlying parent 
lattice is lost. However, within the CALPHAD method- 
ology, one has to define the Gibbs energy of a phase 
throughout the entire composition range, regardless of 
whether the structure is stable or not. In these cases, 
it is necessary to constrain the relaxations so that they 
are consistent with the lattice vectors and atom positions 
of an hep lattice. Obviously, the energetic contributions 
due to local relaxations are not considered in this case. 
The results of these constrained relaxations can there- 
fore be directly compared to those calculations using the 
CPA. In most cases, local relaxations were not significant. 
However, in a few instances, it was found that the struc- 
ture was too distorted to be considered as hep after the 
full relaxation. However, this symmetry check was not 
sufficient to characterize the relaxation behavior of the 
relaxed SQS. Furthermore, in some of the cases it may 
be possible for the structure to fail the symmetry test 
and still retain an hcp-like environment within the first 
couple of coordination shells, implying that the energet- 
ics and other properties calculated from these structures 
could be characterized as reasonable, although not opti- 
mal, approximations of random configurations. 



1. Radial Distribution Analysis 

In order to investigate the local relaxation of the fully 
relaxed SQS, their radial distribution (RD) was analyzed. 
Through this analysis, the bond distribution and coordi- 
nation shells were studied to determine whether the re- 
laxed structures maintained the local hcp-like environ- 
ment they were supposed to mimic in the first place. 
Additionally, this analysis permitted us to quantify the 
degree of local relaxations up to the fifth coordination 
shells. 

The RD of each of the fully relaxed structures was ob- 
tained by counting the number of atoms within bins of 
10~'^A, up to the fifth coordination shell. In order to 
eliminate high frequency noise, the raw data was scaled 
and smoothed through Gaussian smearing with a charac- 
teristic distance of 0.01 A. Pseudo-Voigt functions were 
then used to fit each of the smoothed peaks and the good- 
ness of fit was in part determined through the summation 
of the total areas of the peaks and comparing them to the 
total number of atoms that were expected within the an- 
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alyzed coordination shells. The relaxation of the atoms 
at each coordination shell is quantified by the width of 
the corresponding peak in the fitted RD. 

The RD results of selected SQS's are given in FIG. H 
The unrelaxed, fully relaxed, and non-locally relaxed 
structures are compared in each case as well as the 
smoothed bond distributions and their fitted curves. 
These results are representative of the RD's obtained for 
the seven binary systems at the three compositions stud- 
ied. 



FIG. |2(a)| shows the RDs for the Hf-Zr SQS at the 50 
at.% composition. As can be seen in the figure, the RDs 
for the unrelaxed and non-locally relaxed SQS are almost 
identical, implying that in this system Vegard's Law is 
closely followed. Furthermore, the RD for the fully re- 



laxed SQS in FIG. 2(b) shows a rather narrow distribu- 



tion around each of the the bond-lengths corresponding 
to the ideal or unrelaxed structure. The system there- 
fore needs to undergo very negligible local relaxations in 
order to minimize its energy. 

In the case of the Cd-Mg solution at 50 at.% 
(FIG. 2(c) I, the RDs of the unrelaxed and non-locally 



relaxed SQS are more dissimilar. Even in the non-locally 
relaxed calculation, the original first coordination shell 
(corresponding to the six first-nearest neighbors) has 
split into two different shells (of 4 and 2 atoms) and the 
position of the peak is noticeably shifted. The first two 
well defined coordination shells of the unrelaxed structure 
have merged into a single, broad peak at 3.14A upon full 
relaxation, as shown in FIG. 2(d) This peak now en- 



closes 12 first nearest neighbors. As shown in TABLE [III 
112,1 and 1124 have two different types of pairs. However, 
since they have the same correlation functions, they can- 
not be distinguished. In FIG. 2(d)| it is also shown how 
the fourth and fifth coordination shells merge at 
enclosing 18 atoms. It can be expected that if the c/a 
ratio of a relaxed structure is close to ideal and the broad- 
ening of nearby shells are wide enough that they merge, 
then the structure has almost the same radial distribu- 
tion of an ideal hep structure, albeit with a larger peak 
width. 



FIG. 2(e) shows the RD for the MgsoZrso composi- 
tion. Among the three RD's presented in FIG. [2l this 
one is clearly the one that undergoes the greatest dis- 
tortion upon full relaxation. Even in the non-locally re- 
laxed structures there is a broad bond- length distribution 
around the peaks of the unrelaxed SQS. With respect to 
the fully relaxed SQS, it can be seen how the peaks for 
the fifth and sixth coordination shells have practically 
merged. In this case, the local environment of each atom 
within the SQS stops being hcp-like within the first cou- 
ple of coordination shells. Although the two end mem- 
bers of this binary alloy have an hep as the stable struc- 
ture, it is evident from this figure that the SQS arrange- 
ment is unstable and there is a tendency for the structure 
to distort. In this system, there is a miscibility gap in 
the hep phase up to ~ 900K and the RD refiects the 
tendency for the system to phase-separate. 



The results from the peak fitting for all the fully re- 
laxed SQS's are summarized in TABLE IIVI It should 
be noted that regardless of the system and compositions, 
the sum of the areas under each peak should converge to 
a single value, proportional to 50 atoms. For each peak, 
the error was quantified as the absolute and normalized 
difference between the expected and actual areas. The 
error reported in the table is the averaged value for all 
the peaks in the RD. The broadness of the peaks in the 
RD is quantified through the full width at half maxi- 
mum, FWHM. In the table, the reported FWHM corre- 
sponds to the average FWHM observed for the coordina- 
tion shells enclosing a total of 50 atoms. Note that the 
alloys with the smallest FWHM are Hf-Zr and Cd-Mg. 
As will be seen later, Hf-Zr behaves almost ideally and 
Cd-Mg is a system with rather strong attractive interac- 
tions between unlike atoms that forms ordered hexagonal 
structures at the 25 and 75 at.% compositions. 



2. Bond Length Analysis 

In addition to the RD analysis, we performed the bond 
length analysis (A-A, B-B, and A-B) for all the relaxed 
SQS's. In TABLE El the bond lengths corresponding to 
the first nearest neighbors for all the 21 SQS's are pre- 
sented. As expected, in the majority of the cases the 
sequence da < dij < djj is observed throughout the 
composition range, where dij corresponds to the bond 
distance between two different atom types. The two no- 
table exceptions to this trend correspond to the Cd-Mg 
and Mg-Zr alloys. As will be mentioned below, the Cd- 
Mg system tends to form rather stable intermetallic com- 
pounds at the 25, 50 and 75 at.% compositions, including 
two hexagonal intermetallic compounds. The calculated 

-is 



enthalpy of mixing in this case — shown in FIG. 3(a) 



the most negative among seven binaries studied and the 
fact that the Cd-Mg bonds are shorter than Cd-Cd and 
Mg-Mg seems to reflect the tendency of this system to 
order. In the case of the Mg-Zr alloys, the Mg-Zr bonds 
are longer than Mg-Mg and Zr-Zr, suggesting that this 
system has a great tendency to phase separate, as indi- 
cated by the presence of a large hep miscibility gap in 
the Mg-Zr phase diagrami^. 



B. Enthalpy of Mixing 

It is obvious that if an hep SQS alloy is not stable 
with respect to local relaxations, its properties are not ac- 
cessible through experimental measurements. However, 
approximate effective properties could still be estimated 
through CALPHAD modeling. In order to compare the 
energetics and properties of the calculated SQS's with the 
available experiments or previous thermodynamic mod- 
els, only the non-locally relaxed structures were consid- 
ered whenever the SQS was identified as unstable. This 
effectively assumes that the structures in question are 
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FIG. 2: Radial distribution analysis of selected SQS's. The dotted lines under the smoothed and fitted curves are the error 
between the two curves. 



constrained to maintain their symmetry. The total ener- 
gies of the structures under symmetry-preserving relax- 
ations are obviously higher since the relaxation energy 
is not considered. However, we can consider these cal- 
culated thermochemical properties as an upper bound 
which can still be of great use when attempting to gen- 



erate thermodynamically consistent models based on the 

combined first-principles/CALPHAD approach. 

As mentioned earlier, obtaining thermodynamic prop- 
erties of random alloys using cluster expansion or the 
CPA method has some drawbacks. These methods, how- 
ever, have the advantage of calculating the properties of 
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TABLE IV: Results of radial distribution analysis for the seven binaries studied in this work. FWHM shows the averaged full 
width at half maximum and is given in A. Errors indicate the difference in the number of atoms calculated through the sum of 
peak areas and those expected in each coordination shell. 



Compositions 




Cd-Mg 


Mg-Zr 


AI-Mg 


Mo-Ru 


Hf-Ti 


Hf-Zr 


Ti-Zr 




FWHM 


0.06 ±0.01 


0.09 ± 0.03 


0.08 ± 0.02 


N/A" 


0.11 ±0.03 


0.02 ± 0.00 


0.16 ±0.05 


A75B25 


Prrnr % 


0.72 


0.39 


0.47 


N /A 


1.07 


1.84 


1.27 




Symmetry 


PASS 


PASS 


PASS 


FAIL 


PASS 


PASS 


FAIL 




FWHM 


0.07 ±0.02 


0.15 ±0.02 


0.15 ± 0.07 


0.13 ±0.01 


0.16 ±0.02 


0.03 ±0.01 


0.19 ±0.06 


A50B50 


Error, % 


0.30 


1.42 


1.28 


1.90 


0.35 


1.84 


2.39 




Symmetry 


PASS 


FAIL 


FAIL 


PASS 


PASS 


PASS 


PASS 




FWHM 


0.04 ± 0.01 


0.09 ±0.03 


0.10 ±0.02 


0.07 ± 0.02 


0.11 ±0.06 


0.03 ± 0.00 


0.13 ±0.07 


A25B75 


Error, % 


2.05 


1.22 


0.26 


1.93 


0.26 


1.01 


0.96 




Symmetry 


PASS 


PASS 


PASS 


PASS 


PASS 


PASS 


PASS 



"The radial distribution analysis of Mo 75 at.%-Ru 25 at.% was 
not possible since it completely lost its symmetry as hep. 



TABLE V: First nearest-neighbors average bond lengths for the fully relaxed hep SQS of the seven binaries studied in this 
work. Uncertainty corresponds to the standard deviation of the bond length distributions. 



Compositions 


Bonds 


Cd-Mg 


Mg-Zr 


Al-Mg 


Mo-Ru 


Hf-Ti 


Hf-Zr 


Ti-Zr 


AiooBo 


A-A 


3.07 


3.18 


2.87 


2.75 


3.13 


3.13 


2.87 




A-A 


3.17 ±0.10 


3.18 ± 0.03 


2.92 ±0.03 




3.14 ±0.05 


3.18 ± 0.03 


2.96 ±0.07 


A75B25 


A-B 


3.16 ±0.11 


3.18 ± 0.05 


2.95 ±0.03 


N/A 


3.10 ±0.05 


3.18 ± 0.03 


3.02 ±0.07 




B-B 


3.18 ±0.10 


3.12 ± 0.10 


2.96 ±0.03 




3.09 ± 0.06 


3.18 ± 0.04 


3.04 ± 0.06 




A-A 


3.16 ±0.04 


3.16 ±0.04 


2.98 ±0.06 


2.81 ± 0.08 


3.09 ± 0.06 


3.18 ± 0.03 


3.00 ± 0.09 


A50B50 


A-B 


3.12 ±0.04 


3.20 ±0.06 


3.02 ±0.06 


2.75 ± 0.04 


3.05 ±0.07 


3.19 ± 0.03 


3.06 ± 0.08 




B-B 


3.15 ±0.03 


3.14 ±0.08 


3.07 ±0.08 


2.75 ± 0.04 


3.00 ± 0.06 


3.20 ±0.03 


3.12 ±0.08 




A-A 


3.16 ±0.01 


3.15 ± 0.04 


3.06 ±0.04 


2.73 ± 0.04 


3.02 ± 0.05 


3.19 ± 0.03 


3.09 ± 0.08 


A25B75 


A-B 


3.14 ±0.02 


3.19 ± 0.04 


3.08 ± 0.04 


2.73 ± 0.04 


3.00 ± 0.06 


3.19 ± 0.03 


3.11 ±0.06 




B-B 


3.15 ±0.01 


3.18 ±0.04 


3.11 ±0.03 


2.71 ±0.04 


2.95 ±0.05 


3.20 ±0.04 


3.17 ±0.06 


AqBioo 


B-B 


3.18 


3.19 


3.18 


2.68 


2.87 


3.19 


3.19 



random alloys at arbitrary and closely spaced concen- 
trations. SQS's in this case are at a disadvantage since 
the size of the SQS itself limits the concentrations with 
random-like correlations. Nevertheless, if we can acquire 
the properties at these three compositions, we can suf- 
ficiently describe the tendency of the system. Further- 
more, these SQS's can be applied directly to other binary 
systems without any modifications. 

The enthalpies of mixing for these alloys were calcu- 
lated at the 25, 50, and 75 at.% concentrations through 
the expression: 



Ai/(Ai_,S,) = E{A^^A)-{l-x)E{A)-xE{B) (2) 

where E{A) and E{B) are the reference energies of the 
pure components in their hep ground state. 

In the following sections, the generated SQS's are 
tested by calculating the crystallographic, thermody- 
namic and electronic properties of hep random solutions 
in seven binary systems, Cd-Mg, Mg-Zr, Al-Mg, Mo-Ru, 
Hf-Ti. Hf-Zr, and Ti-Zr. The results of the calculations 



are then compared with existing experimental informa- 
tion as well as previous calculations. 



C. Cd-Mg 

In the Cd-Mg system, both elements have the same 
valence and almost the same atomic volumes. Conse- 
quently, there is a wide hep solid solution range as well 
as order/disorder transitions in the central, low temper- 
ature region of the phase diagram. In fact, at the 25 
and 75 at.% compositions there are ordered intermetallic 
phases with hexagonal symmetries. 

compares the enthalpy of mixing calculated 
from the fully relaxed and symmetry preserved SQS with 
the results from cluster expansion^. The results by Asta 
et al.-^ at 900K are presented for comparison since it is 
to be expected that these values would be rather close to 
the calculated enthalpy of completely disordered struc- 
tures. The previous and current calculations are also 
compared with the experimental measurements as re- 
ported in Hultgren-- at 543K. The first thing to note 
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from FIG. 3(a) is that the fully relaxed and symmetry 
preserved calculations are very close in energy, implying 
negligible local relaxation. Additionally, the present cal- 
culations are remarkably close ('^1 kJ/mol) to the exper- 
imental measurements. By comparing the SQS enthalpy 
of mixing with the results from the cluster expansion cal- 
culations^, it is obvious that the former is, at least in 
this case, more capable of reproducing the experimental 
measurements . 

Formation enthalpies of the three ordered phases in 
the Cd-Mg system, CdgMg, CdMg, and CdMgg are also 
presented. The measurements from Hultgren— deviate 
from the calculated results from Asta et al.— and this 
work. Cd and Mg are known as very active elements 
and it is likely that reaction with oxygen present during 
the measurements may have introduced some systematic 
errors. Furthermore, the measurements were conducted 
at relatively low temperatures, making it difficult for the 
systems to equilibrate. Nevertheless, experiments and 
calculations agree that these three compounds constitute 
the ground state of the Cd-Mg system. 



FIG. 3(b) also shows that the present calculations are 



able to reproduce the available measurements on the vari- 
ation of the lattice parameters of hep Cd-Mg alloys with 
composition, as well as the deviation of these parameters 
from Vegard's Law. This deviation is mainly related to 
the rather large difference in c/a ratio between Cd and 
Mg. The c/a ratio of Cd is one of the largest ones of all 
the stable hep structures in the periodic table. 



D. Mg-Zr 

The Mg-Zr system is important due to the grain refin- 
ing effects of Zr in magnesium alloys. According to the 
assessment of the available experimental data by Nayeb- 
Hashemi and Clark—, the Mg-Zr system shows very little 
solubility in the three solution phases, bcc, hep and liq- 
uid. In fact, the low temperature hep phase exhibits a 
broad miscibility gap up to 923K, corresponding to the 
peritectic reaction hep + liquid — > hcp^. 

Our calculations yielded a positive enthalpy of mixing, 
confirming the trends derived from the thermodynamic 
model developed by Hamalainen et al.— . In the case 
of the full relaxation, however, it was observed that the 
Mg5oZr5o SQS was unstable with respect to local relax- 
ations. The instability at this composition and the large, 
positive enthalpy of mixing indicate that the system has 
a strong tendency to phase-separate. By comparing the 
fully relaxed and the non-locally relaxed structures, we 
estimate that the local relaxation energy lowers the mix- 
ing enthalpy of the random hep SQS by about 2 kJ/mol 
in this system. 

FIG. [3] shows the calculated mixing enthalpy for the 
Mg-Zr hep SQS with no local relaxations, as well as 
the mixing enthalpy calculated from the thermodynamic 
model by Hamalainen et al.-2&, which was fitted only 
through phase diagram data. It is therefore remark- 



able that the maximum difference between the CAL- 
PHAD model and the present hep SQS calculations is 
^3 kJ/mol. The CALPHAD model, however, docs not 
correctly describe the asymmetry of the mixing enthalpy 
indicated by the first-principles calculations. The results 
of the hep SQS calculations for the Mg-Zr system have 
recently been used to obtain a better thermodynamic de- 
scription of the Mg-Zr system^ and, as can be seen in the 
figure, this description is better at describing the trends 
in the calculated enthalpy of mixing. 



E. Al-Mg 



As one of the most important industrial alloys, the Al- 
Mg system has been studied extensively recently22i2^i2^. 
This system has two eutectic reactions and shows sol- 
ubility within both the fee and hep phases. However, 
the solubility ranges are not wide enough so there is only 
limited experimental information for the properties of the 
hep phase. The maximum equilibrium solubility of Al in 
the Mg-rich hep phase is around 12 at.%. 



In FIG. 5(a) the calculated enthalpy of mixing is 



slightly positive. The fully relaxed calculations show that 
the SQS with the 50 at. % composition was unstable with 
respect to local relaxations. This can be explained by the 
strong interaction between Al and Mg, as evident from 
the tendency of this system to form intermetallic com- 
pounds at the middle of the phase diagram, such as /3- 
Ali4oMg89, 7-Ali2Mgi7, and e-AlgoMgas. At the 25 and 
75 at.% compositions the SQS's were stable with respect 
to local relaxations because both elements have a close- 
packed structure. Furthermore, at these compositions ei- 
ther the fee or hep phase take part in equilibria with some 
other (intermetallic) phase. FIG. 5(a) shows that the 



present fully relaxed calculations are in excellent agree- 
ment with the most recent CALPHAD assessments^LS^. 
Note also that in this case, and contrary to what is ob- 
served in the Cd-Mg binary, the energy change associated 
with local relaxation is not negligible, although it is still 
within ~1 kj/mol. 

Additionally, the calculated lattice parameters agree 
very well with the experimental measurements of Mg- 
rich hep alloys, as can be seen in FIG. 5(b) It is im- 



portant to note that the lattice parameter measurements 
of metastable hep alloys from Luo et al."22.(77.4 and 87.8 
Mg at.%) are lying on the extrapolated line between the 
75 at.% SQS and the pure Mg calculations. This is an- 
other example of how SQS's can be successfully used in 
calculating the properties of an hep solid solution system 
with narrow solubility range and mixed with non-hep el- 
ements, even in the metastable regions of the phase dia- 
gram. 
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FIG. 3: Calculated and experimental results of mixing enthalpy and lattice parameters for the Cd-Mg system 



14 



12 



SQS:Fully Relaxed 
SQS:SyiTuiietry Preserved 
1999 Hamalainen et al.(assessment) 
2005 Arroyave et al.(assessment) 




0.25 0.5 0.75 

Mole Fraction, Zr 



In the Ru-rich side, the maximum solubility of Mo in 
the hcp-Ru matrix is up to 50 at.%. The calculations at 
M025RU75 and M050RU50 retained the original hep sym- 
metry but M075RU25 did not. The instability of the Mo- 
rich SQS is not surprising since the Mo-rich bee region 
is stable over a wide region of the phase diagram. As 
shown in Wang et al.— , elements whose ground state is 
bcc are not stable in an hep lattice and viceversa (bee 
Ti, Zr and Hf are only stabilized at high temperature 
due to anharmonic effects) . Thus hep compositions close 
to the bcc-side would be dynamically unstable and would 
have a very large driving force to decrease their energy 
by transforming to bcc. 



FIG. 4: Calculated enthalpy of mixing in the Mg-Zr sys- 
tem compared with a previous thermodynamic assessmenls^. 
Both reference states are the hep structure. 



F. Mo-Ru 

The Mo-Ru system shows a wide solubility range 
within both the bcc and hep sides of the phase diagram. 



Recently, Kissavos et al. - calculated the enthalpy of 
mixing for disordered hep Mo-Ru alloys through the CPA 
in which relaxation energies were estimated by locally re- 
laxing selected multi-site atomic arrangements. Enthalpy 
of formation for hep solutions were calculated from Eqn[2] 
shown below. The enthalpy of mixing of the disordered 
hep phase can be evaluated accordingly based on the so- 
called lattice stabihty2, E'""'{Mo) - E'"'p{Mo). 
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FIG. 5: Calculated and experimental results of mixing enthalpy and lattice parameters for the Al-Mg system 



- E'''P{Mo,-^Ru,) - (1 
-(1 - x)E'""'{Mo) + (1 



x)E^'"'{Mo) - xE'"'P{Ru) 
x)E'"'P{Mo) - xE'"'P{Ru) 
- x)E'"'P{Mo) 

x)[E^'"'{Mo) - E^-'PiMo)] (3) 



Usually, structural energy differences (or lattice stabil- 
ity) between first-principles calculations and CALPHAD 
show quite good agreement. However, for some transition 
elements, the disagreement between the two approaches 
is quite significant^^. Mo is one such case, with the struc- 
tural energy difference between bcc and hep from first- 
principles calculations and the CALPHAD approach dif- 
fering by over 30 kJ/mol. After a rather extensive anal- 
ysis, Kissavos et al. arrived at the conclusion that in 
order to reproduce enthalpy values close enough to the 
available experimental data'^^ the CALPHAD lattice sta- 
bility (11.55 kJ/mol) needed to be used for the value of 
the bcc hep promotion energy. 

The SQS and CPA calculations are compared with the 
experimental measurements in FIG. [51 On the assump- 
tion that the experimental measurements by Kleykamp^ 
are correct, the derived enthalpy of formation of the hep 
Mo-Ru system from the first-principles calculated lattice 
stability with the SQS and CPA approach in FIG. 6(a) 
cannot reproduce the experimental observation at all 
since the first-principles bcc hep lattice stability for 
Mo is 42 kJ/mol. Given this lattice stability, the only 
way in which the first-principles calculations within both 



the SQS and CPA approaches would match the exper- 
imental results would be for the calculated enthalpy of 
mixing to be very negative, which is not the case. In 
fact, as can be seen in FIG. 6(a) the SQS and CPA cal- 



culations are very close to each other. 



On the other hand, the enthalpy of formation derived 
from the CALPHAD lattice stability in FIG. |6(b)| shows 
a better agreement than that from the first-principles 
lattice stability. It is important to note that the CAL- 
PHAD lattice stability was obtained through the extrap- 
olation of phase boundaries in phase diagrams with Mo 
and stable hep elements and, therefore, are empirical. 
The reason why such an empirical approach would yield 
a much better agreement with experimental data is still 
the source of intense debate within the CALPHAD com- 
munity and has not been resolved as of now. The main 
conclusion of this section, however, is that the SQS's were 
able to reproduce the thermodynamic properties of hep 
alloys as good as or better than the CPA method while 
at the same time allowing for the ion positions to locally 
relax around their equilibrium positions. 
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FIG. 6: Enthalpy of formation of the Mo-Ru system with both first principles and CALPHAD lattice stabilities, 
states are bcc for Mo and hep for Ru. 



Reference 



G. IVA Transition Metal Alloys 

The group IVA transition metals, Ti, Zr, and Hf have 
hep structure at low temperatures and transform to bcc 
at higher temperatures due to the effects of anharmonic 
vibrations. When they form a binary system with each 
other, they show complete solubility for both the hep 
and bcc solutions without forming any intermetallic com- 
pound phases in the middle. 

The Hf-Ti binary is reported to have a low temperature 
miscibility gap and was modeled with a positive enthalpy 
of mixing by Bittermann and Rogl"^^. FIG. |7(a) shows 
remarkable agreement between the fully relaxed first- 
principles calculations and the thermodynamic model, 
which was obtained by fitting the experimental phase 
boundary data. Despite the fact that the local relax- 
ation energies are rather large (~ 4 kJ/mol), the lattice 
parameters in both cases agree between each other and 
with the experimental results^Sii^i^. 

In the case of the Ti-Zr binary, although no low- 
temperature miscibility gap has been reported, Kumar 
et al.— found that the enthalpy of mixing for the hep 
solutions in this binary was positive through fitting of 
phase diagram data. Our results confirm this finding, 
although with even more positive enthalpy. They are in 
fact similar in value to those calculated in the Hf-Ti al- 
loys, suggesting that a low temperature miscibility gap 
may also be present in this binary. 

In the Hf-Zr system no miscibility gap has been re- 
ported. The hep phase was modeled as an ideal solu- 
tion {AHmix = 0) in the CALPHAD assessmenl^^i. The 



present calculations suggest that the enthalpy of mixing 
of this system is positive, although rather small. In this 
case, it is expected that any miscibility gap would only 
occur at very low temperatures. 

The three systems described in this section are chem- 
ically very similar, having the same number of electrons 
in the d bands. Electronic effects due to changes in the 
widths and shapes of the DOS of the d bands are not 
expected to be significant in determining the alloying en- 
ergetics. Charge transfer effects are also expected to be 
negligible. The enthalpy observed can then be explained 
by just considering the atomic size mismatch between the 
dilferent elements. As was shown in TABLE |Vl the Hf- 
Zr hep alloys are the ones with the smallest difference in 
their lattice parameter, thus explaining their very small 
positive enthalpy of mixing. 

As a final analysis of the ability of the generated SQS 
to reproduce the properties of random hep alloys, FIG.[H] 
shows the alloying effects on the electronic DOS in Ti-Zr 
hep alloys. The figure also presents the results obtained 
through the CPA approach by Kudrnovsky et al.— . As 
can be seen in the figure, both calculations predict that 
the DOS corresponding to the occupied d states are vir- 
tually insensitive to alloying. The overall shape of the d- 
DOS remains relatively invariant. Since Ti and Zr have 
the same number of valence electrons, the fermi level re- 
mains essentially unchanged as the concentration varies 
from pure Zr to pure Ti. On the other hand, alloying 
effects are more pronounced in the c?-DOS corresponding 
to the unoccupied states. FIG. [5] shows how the broad 
peak at ^ 4.5 eV of the d-DOS for Zr is gradually trans- 
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FIG. 7: Enthalpy of mixing for the Hf-Ti, Hf-Zr and Ti-Zr binary hep solutions calculated from first-principles calculations and 
CALPHAD thermodynamic models. All the reference states are hep structures. 



formed into a narrow peak at ~ 3.0 eV as the Ti content 
in the alloy is increased. The results from the CPA and 
the first-principles SQS calculations thus agree with each 
other, confirming the present results. 



V. SUMMARY 

We have created periodic special quasirandom struc- 
tures with 16 atoms for binary hep substitutional alloys 
at three dilferent compositions, 25, 50, and 75 at.%, to 



mimic the pair and multi-site correlations of random so- 
lutions. 

The generated SQS's were tested in seven different bi- 
naries and showed fairly good agreement with existing ex- 
perimental either enthalpy of mixing and/or CALPHAD 
assessments and lattice parameters. Analysis of the ra- 
dial distribution and bond lengths in the 21 calculated 
SQS's, yielded a detailed account of the local relaxations 
in the hep solutions and has been proven a useful way of 
characterizing the degree relaxation over several coordi- 
nation shells. 
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FIG. 8: Calculated DOS of Tii^^Zrx hep solid solutions from (a) SQS and (b) CPA^ 



It should also be noted that when using enthalpy of 
mixing to derive formation enthalpy to compare with ex- 
perimental measurements, there can be a severe discrep- 
ancy between theoretical calculations and experimental 
data when the lattice stability, or structural energy dif- 
ference, from first-principles calculation is problematic 
such as the Mo-Ru system in this work. This problem 
remains as an unsolved issue. 

These supercells can be applied directly to any substi- 
tutional binary alloys to investigate the mixing behavior 
of random hep solutions via first-principles calculations 
without creating new potentials, as in the coherent po- 
tential approximation (CPA) or calculating other struc- 
tures in the cluster expansion. Although the size of the 
current SQS's is not large enough to generate a super- 
cell which can satisfy its correlation function at more 
than just three compositions (x = 0.25, 0.5, and 0.75 in 



binary), calculations for these compositions can 
yield valuable information about the overall behavior of 
the alloys. 
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